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Summary 

Lifetimes of low-altitude lunar orbits are studied 
in this report to identify feasible parking orbits for 
future lunar missions. The lunar missions currently 
under study, unlike the Apollo missions, involve long 
stay times. To determine orbital lifetimes of lunar 
parking orbits, a model to describe the nonspherical 
mass distribution of the Moon must be adopted. A 
short discussion of prfnhous attempts to describe the 
Mooirs gravitational field is included, and emphasis 
is placed on spherical harmonic-gravity models. A 
subset of the hfth-order and fifth-degree lunar gravity 
model was adopted in this investigation to generate 
orbital lifetime predictions. This simplified model 
consists of the five gravitational coefficients J 2 , 

<^ 22 * C:i\. 

The primary perturbation on a low-altitude 
(100 km or 300 km) lunar parking orbit is the Moorrs 
nonspherical gravity fi('ld. In this analysis, all other 
perturbations (third body perturbations, solar radi- 
ation i)ressure, etc.) ar(' iK'glected. Although per- 
turbations due to the Earth and the Sun are not 
negligible for lunar orbits, theses effects are small com- 
pared with tfie Moon's nonspherical potential for all 
the orbits considered in this analysis. By investi- 
gating the effects that initial conditions have on the 
subsequent lifetime of an orbit, a technique is intro- 
duced to aid mission planners in the selection of lunar 
parking orbits. In this investigation, tlie lifetimes of 
near-circular parking orbits, with various initial or- 
bital elements and either l()()-km or 300-krn perilune 
altitude, are analyzed for mission planning purposes. 

Orbital lifetimes are heavily dependent on the 
initial conditions of the orbit, particularly the ini- 
tial inclination and argument of perilune. The lu- 
nar gravity modc'l utilized in this analysis for the 
100-km initial perilune altitude case yields lifetime 
predictions of less tfian 40 days for some orf)its, and 
more than a year for others. Five distinct bands of 
short-lifetimc' orbits a])pear as a function of tlu' ini- 
tial inclination; tlu'se bands are separated by bands 
of long-lifetime orlhts. Of particular interest is a set 
of orl)its with an inclination of approximately 70°; 
this set of orbits yields long lifetimes and provides 
the high latitude coverage that is desirable for vari- 
ous missions. Tlie »/5 coefficient contributes the? dom- 
inant effect in perilune altitude decay and, therefore, 
orbital lifetiim's. 

The methods presenU'd in this report are suit- 
able for incorporating the Moon s nonspherical grav- 
itational effects into the preliminary design level for 
future lunar mission planning. However, inconsisten- 
cies and limitations, caused primarily by a lack of 


satellite tracking data from the far side of the Moon, 
are inherent in all existing lunar gravity models. The 
uncertainty in orbital lifetime predictions due to er- 
rors in the lunar gravity model is addressed through 
the use of sensitivity coefficients. The uncertainty in 
the rate of perilune altitude decay that corresponds 
to the uncertainty in the values of the coefficients 
for the gravity model adopted in this analysis is pre- 
sented. Also, plots of the values of the sensitivity 
coefficients, which can be used to evaluate the uncer- 
tainty in perilune-altitude decay rates of each grav- 
itational coefficient for any lunar gravity model, are 
presented. 

Introduction 

President Bush's proposal of a Space Exploration 
Initiative in 1989 sparked a renewed interest in lunar 
mission planning. The objectives outlined in this ini- 
tiative include the establishment of a permanent lu- 
nar outpost. (See refs. 1 and 2.) Lunar stay times on 
the order of 30 to 180 days will be required for initial 
deployment and ongoing support of the outpost. In 
some analyses, preliminary missions will involve the 
insertion of a satellite into a low-altitude lunar orbit 
to map the Moon’s surface. Also in these analyses, 
initial manned missions will require tlu^ placement 
of a lunar transfer vehicle in a low-altitude parking 
orbit. Ongoing support of an outpost might include 
the placement of a space station in low lunar orbit. 
Any of these missions will requiri' spacecraft to be in 
lunar parking orbits for long periods of time. Pre- 
vious orbital determination studies for lunar satel- 
lites have indicated that the Moon’s nonspherical 
gravity field will have a great effect on the sub- 
sequent lifetime of the orbit. (See refs. 3 to 5.) In 
fact, a lunar-orbit space station i)roposed by NASA 
to be in a 60-n.mi. (Ill km) circular polar orbit 
about the Moon was found to impact the lunar sur- 
face in 140 days if no station-keeping altitude boost 
maneuvers were performed. (See ref. 6.) 

An important element involved in establishing a 
base on the Moon is the initial manned mission to 
the lunar surface. As with the Apollo missions, the 
lunar transfer vehicle (LTV) will be inserted into a 
parking orbit about the Moon at arrival. Apollo mis- 
sions restricted these parking orbits to circular, near- 
equatorial orbits. However, the utilization of park- 
ing orbits at all inclinations must be addressed to 
accommodate the vast array of lunar missions cur- 
rently being proposed. The LTV will remain in this 
parking orbit until a departure burn for return to 
Earth is initiated. The transh'r of cargo and person- 
nel to the lunar surface is accomplished by the de- 
scent of a lunar lander from the parking orbit. Ascent 
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Figui'(' 1. Eartli- lunar mission sequenco. (From ref. 1.) 

and rendezvous of the lander with the LTV in the 
parking orbit initiates the Earth-return sc^quence of 
the mission. (See fig. 1.) Therefore, long stay times 
on the lunar surface will require that the LTV remain 
in a parking orbit for a long period. The mission 
designer must address the variations in the orienta- 
tion of the parking orbit for the duration of the mis- 
sion to establish tlu' necessary AV requir(unents for 
performing the rendezvous and trans-Earth injection 
maneuvers. 

Since the A{)ollo missions involved only short stay 
times at the Moon (on the order of 1 to 3 days), little 
work w^is perfornuHl on predicting long term changes 
in the parking orbits of lunar modules. However, for 
longer stay times (30 to 180 days), an accurate model 
of the Moon's gravity field is required for preliminary 
mission planning, especially with low-altitude park- 
ing orbits. The need for long orbital lifetimes is a key 
criterion in the process of selecting feasible parking 
orbits. By investigating the effects that initial con- 
ditions have on the subsequent lifetime of an orbit, a 
technique is introduced that will aid mission planners 
in the selection of lunar parking orbits. 

Be'cause of the Moon’s irregular internal structure 
and surface shape, the lifetime of a lunar satellite or- 
bit is constrained by the nonspherical nature of the 
Moon's gravity field. In the present analysis, a sim- 
plified gravitational model of the Moon is introduced 
that will enable mission planners to easily predict 
long-term changes in lunar parking orbits at the pre- 
liminary design level. The development of a simpli- 
fied gravitational model with stifficient accuracy is 
necessary to invt'stigate orbital lifetimes for the large 
number of initial orbital pararnetc'rs possible. Uti- 
lization of a simplified model wnll significantly reduce 
tlu' reciuired (‘omputational time needed to perform 
this analysis. The effect of individual terms in grav- 
ity mod(‘ls on orbital lifetimes is also addressed. By 
investigating the effects of the lunar gravity model on 
tlu' various parking orbits, the parameters that are 
most important in determining lifetime predictions 
are identified. 
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Fourier coefficients 
semimajor axis, km 
gravitational coefficients 
infinitesimal mass clement, kg 
eccentricity 

eccentricity averaged over period of 
third body about central body 

new value of eccentricity (after 
integration) 

old value of eccentricity (before 
integration) 

gravitational constant, 

6.67 X 10 m'*-kg 

altitude, km 
inclination, deg 

inclination averaged over period of 
third body about central body, deg 

zonal gravitational coefficient, -C„o 

lunar transfer vehicle 

mean anomaly, deg 

mass of Moon, 7.35 x 10^^ kg 

mean angular motion, rad/day 

mean angular motion averaged over 
period of third body about central 
body, rad/day 

mean angular motion of third body 
about central body, rad/day 

associated Legendre polynomial 

position vector 

radius of Moon, 1739 km 

distance between mass (dement and 
exterior point, km 

radius, km 

rate of altitude decay, km/day 
time, days 

integration time step, 1 day 
gravitational potential 
gradient of gravitational potential 
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V" 

disturbing function (nonspherical 
gravitational potential) 

AV 

velocity increment, m/sec 

V 

volume, km^ 

w 

rotation rate of Moon, 13-18°/day 

J/, 2 

position in Cartesian coordinates, km 

i,y,z 

velocity in Cartesian coordinates, 
km / sec 

A 

longitude, deg 

1/ 

true anomaly, deg 

P 

density, g/enV^ 

nm 

standard deviation of gravitational 
coefficient 

(Jf 

0 ' 

variance in perilune altitude decay 
rate, km/day 

<t> 

latitude, deg 

n 

longitude of ascending node, deg 

iii 

inertial node longitude, deg 


selenographic node longitude, deg 


^ contribution for C 31 term 

ijj 

argument of perilune, deg 

UJ 

argument of perilune averaged over 
period of third body about central 
body, deg 

Subscripts: 


m 

order of gravitational coefficient 

n 

degree of gravitational coefficient 

P 

perilune 

Background 


Gravitational Potential Theory and 
Resulting Perturbations 

To determine the orbital lifetime of a satellite, the 
gravitational field of the attracting body must first 
be described. Since gravity is a conservative force, 
the gravitational field of a body can be represented 
by a potential function. A solution for the form of 
the potential of a body can be obtained in terms of 
a mass integral definition or by finding a solution 
to Laplace’s equation. In the first approach (refs. 7 
and 8 ), infinitesimal mass elements are integrated 
over the entire body to describe the potential at 



Figure 2. Solution for potential of a bo<iy using a mass 
integral definition. 

an exterior point, p(r, (/), A), where r is the radial 
distance between the point and the origin, and (p 
and A are tlie latitude and longitude of the point 
(fig. 2 ) as follows: 

f f dm 

U(r.4>.\)= dU = G ^ 

7 Body ./Body 

=g/ (1) 

where IZ is the distance between the mass element 
and the exterior point p, p is the density function 
of the body, and V is the volume of the body. The 
term is expanded in a Legendre polynomial series, 
and the addition theorem for Legendre associated 
polynomials is used to write the potential in the form 

^ m ( yi+T ) C^m(w^ (f>) 

n=() 

X [Anm COS (mA) + Burn (mA)] (2) 

This spherical harmonic expansion of the i)otential 
function can also be obtained by using a separation- 
of- variables method to solve Laplace’s equation. This 
approach is outlined in reference 9. 

The A{)[) term represents the uniform spherical po- 
tential contribution Assuming that the origin of 
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tlu' (‘oordinate systeiu is at the center of mass of the 
body, then /Iiq = /tu = B\i — 0. By also introduc- 
ing the variables C,im = and Sum = the 

potential can be written in the traditional form as 
follows: 


U = 


GM 


1 + EE 


u=2 /n =() 


n 

Pn7n(sill0) 


X cos {tnX) + 5n/;i sill (mA)] 


(3) 


If the density function of the body is known, the value 
of the coeflidents Cjim and Sum can be deterniiiKKl 
by integration of the appropriate nuiss integrals over 
the volume of the body. In praeti(x\ since the 
density function is not known, the values of the Cum 
and Sum coefficients are determined empirically by 
tracking satellites orbiting the body. Using statistical 
nu'thods, the best set of Cum ^nd Sum coefficients 
that desc‘ribe the orbit mv determined. (See section 
‘'Lunar Gravity Models.'") 

Tlu’ disturbing function \\ defined as V = 
U — contains tlu' nonspherical gravitational 

contribution to tlu' potential. This nonspherical 
gravitational contribution arises from the noii- 
spherical mass or density distribution of th(^ body. 
The net effect of this nonuniformity is the creation 
of a small disturbance force on an orbiting body; 
this disturbance causes a change in the orbital ek^- 
ments oxer time. Tlu' Lagrange planetary equations 
(ref. 10) (k'sc'ribe the effects of the disturbing func- 
tion on an orbiting body as follows (the n term in 
eq, ( If) is a result of the term): 
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These equations are the result of a direct coordi- 
nate transformation of the equations of motion from 
Cartesian coordinates (x, y, z, x, y, z) to t fie six clas- 
sical orbital elements (a, c, i, 0, u?, Af). Once the 
nonspherical gravitational contrilmtion of tlu‘ poten- 
tial V is defined, the partial derivatives required for 
the solution of the Lagrange planetary eciuations can 
be evaluated; this evaluation (uiables the orbit of the 
spacecraft to be determirK'd. The resulting equations 
for several of the nonspherical potential terms arc 
contained in appendix A. These first-order analyti- 
cal equations are derived with an averaging technique 
that involves integration of the .short-pexiod effects in 
the disturbing function. (See refs. 3 and 11.) 

Mass Concentrations — An Alternative 

Approach 

The spherical harmonic C!xpression for the poten- 
tial (ecj. (3)) is very useful for dc'seribing gravit ational 
fields that vary only slightly from a spherical field in 
a smooth manner, such as for the Earth. A detailed 
description of tfu' Earth’s potential is obtained with 
only a few harmonic terms. The Moon is smaller and 
less massiv(' than tlu' p]arth; ther('fon\ it can support 
more stress per unit mass and has tlu' capability to 
possess greater gravitational anomali(^s than does the 
Earth. Bt'cause of its irregular intexnal structure and 
surface shape, the Moon has a very com])li(‘at(Mi grav- 
itational field. As a result of the^se propeudies, the' 
Moon s gravitatie:>nal field caniicjt be^ de^scribexi with 
only a few terms. 

Several problems are ence>untere;el when a spheri- 
cal harmonic e^xpansion of the potential is used. One 
problem is the slenv convergeaice of the exi>ansion for 
points ne'ar the lunar surfae'e, when' y is slightly k'ss 
than e>ne'. For this tease) n. many gravitational co- 
effie:ients arc rcquirexl to (l(^scril)e' the' orbit e>f a k)w- 
altituele satellite. Also, spherical harmonic expan- 
sions are unable to ek'seribe locali/nd gravitational 
anomalie^s unless a large numbeu' of ce)e’fficients are' 
introduced; they are more appropriate* for describing 
an average gravitational field. For example, an ex- 
pansion on the orelcr of 180 (rn — 180, se.*e eej, (3)) is 
necessary to describe a local anomaly that subtends 
an angle of 1°. (See ref. 12.) Difficulties in obtain- 
ing values e)f C^u/n and Sum f^^ the* Moe)ii have^ also 
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been encountered for reasons that are discussed in 
the section “Lunar Gravity Models.'' 

The existence of lunar mass concentrations (mas- 
cons) was postulated in reference 13. These are lo- 
calized regions of higher than average density that 
produce measurable gravity anomalies. A gravity 
anomaly is defined as the residual gravity effect after 
the attraction of a reference body (such as a homo- 
geneous spheroid, rotational ellipsoid, triaxial ellip- 
soid, etc.) is subtracted from the measured gravity 
data (ref. 14). This gravity anomaly is a result of 
the mascon's higher density (mascons have a den- 
sity of approximately 3.3 g/cm*^ whereas surround- 
ing rocks have an average density of only 3.0 g/cm'^). 
Mascons constitute between 0.01 and 0.03 percent 
of the Moon's mass, and their location is coincident 
with the location of lunar seas (maria). (See ref. 12.) 
These gravitational anomalies are found primarily 
along the equator on the near side of the Moon. 

Attempts have been made to use niascon mod- 
els to describe the features of orbital tracking data. 
Models that distribute mascons below the lunar sur- 
face rather than on the surface (ref. 15) and mod- 
els that treat mascons as circular disks rather than 
point masses (ref. 16) have had much more success. 
A model was also developed that treated the niascon 
as a third body to investigate the short-period, long- 
period, and secular effects of mascons on the orbit 
of a spacecraft. (See ref. 17.) The main restriction 
of this latter approach is that it is not ai)plicable to 
low-altitude satellites. 

One advantage of a spherical harmonic model for 
lifetime studies is that it can be used to provide a 
relatively simple analytical approach. Since a spher- 
ical harmonic description allows one to average short- 
period effects (perturbations that result from the 
variation of the mean anomaly M around the orbit), 
these averaging effects can be accounted for over the 
period of a few orbits. Averaging effects for niascon 
models can be developed by mapping the mascons 
to a set of spherical harmonic coefficients. However, 
this process effectively eliminates the advantage of 
the niascon approach, because a large number of co- 
efficients arc needed to accurately model the effects 
of mascons. 

An alternative method was develojied by Ananda 
to construct a disturbing potential in terms of var- 
ious parameters of the mascons and the orbital el- 
ements. (See ref. 18.) This potential is averaged 
over an orbit to eliminate the short-period terms. 
Much of this integration cannot be performed analyt- 
ically. The averaged potential is determined by nu- 
merically evaluating these expressions with Gaussian 


quadrature formulas that require considerable com- 
puter time. Once the potential is defined, the average 
rates of the mean orbital elements can be determined 
by evaluating the derivatives in the Lagrange plane- 
tary equations (eqs. (4)). Work is being conducted 
to develop a method that allows for the application 
of mascon models without the sacrifice of computer 
time. (See ref. 19.) 

Hybrid models, with both spherical harmonics 
and mascons, have also been proposed (refs. 16 
and 20). Purely mascon models, with a large number 
of mascons (on the order of 100), require significant 
computational time. Using a few low-degree spher- 
ical harmonic terms to describe global iionspherical 
contributions decreases the number of mascons re- 
quired to give a complete description of the lunar 
gravity field. The function of the mascons in a hy- 
brid model is only to describe localized gravitational 
effects. This may be the most reasonable gravita- 
tional model to use, because it combines the benefits 
of spherical harmonic models and the mascon ap- 
proach. However, the application of a hybrid model 
on the preliminary mission design level may prove to 
be infeasible because of the excessive computer time 
demanded by the additional complexity associated 
with the mascons. 

A spherical harmonic model was adopted for the 
study of orl)ital lifetimes in this investigation because 
of the availability of these models. Few mascon mod- 
els are available, and metho<ls for including the ef- 
fects of these anomalies in long-term orbital predic- 
tions have not been widely developed. Also, since 
lifetime studies involve averaging effects of the en- 
tire gravitational field of the Moon, a spherical har- 
monic model seems to be more appropriate for this 
application. Mascon models are more appropriate 
for situations in which information about localized 
fields is desired. Some applications might inc:hide 
performing a maneuver with a low-altitude satellite 
above a gravitational anomaly, or calculation of de- 
scent and ascent trajectories to and from a landing 
site near a mascon. If localized gravitational effects 
arc found that cause significant changes in the or- 
bital elements (either as single events or integrated 
over time), spherical harmonic representations of the 
lunar gravity field will need to be abandoned in favor 
of mascon models. 

Lunar Gravity Models 

Derivation of models. Spherical harmonic lu- 
nar gravity models contain values of the gravita- 
tional coefficients Cnm and Snm that define the non- 
spherical contributions to the potential field. (Sec 
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Figure 3. Various types of harmonic coefficients. White 
area r('presents elevation above and black area represents 
elevation below a mean spherical surface. (From ref. 21.) 



ecj. (3).) There are three types of harmonic co- 
efficients; these types are defined by the values of 
the subscripts n and rn. For zonal harmonic co- 
efficients, rn = 0 (Cno, often denoted by — Jn )- These 
coefficients describe an axially symmetric potential 
(in this case, the spin axis of the Moon) indepen- 
dent of tlu' longitude at which the potential is mea- 
sured. Gravity models for the Earth consisting solely 
of zonal harmonics yield good approximations to the 
actual field. However, a longitude dependence in the 
Moon’s gravity field motivates the need to include 
other types of coefficients in lunar models. This de- 
pendence' is introduced by the presence of sectorial 
and tesseral harmonic coefficients. Sectorial harmon- 
ics (where n ~ rn ^ 0) give rise' to zero values of the 
pe)tential only along mcrieiians of longitude, wherciiis 
tesseral harmonics (where n ^ m, and m and n are 
nonzere)) also give rise to zero valuers along parallels 
of latitude. (Se'e' fig. 3, from red. 21.) Once the Cum 
and Sfini ('oefficients arc spe^cified, the potential is 
(‘omplete'ly eleterniined as a superposition of the 
inelivielual harmonic terms. 

Further insight into the physical significance of 
the' individual gravitatiemal-coefficient values can be 
e)btained by ititer preding the coefficients as surface 
deviations from a he)mogeneous sphere. The J 2 gravi^ 
tational e*oefficie'nt analytically repre^sents the oblate- 
ness of a boely (this eeiuatenial ''bulge,'' comme:)n 
to all rotating boelies, arises from the “centrifugal 
force” prodiu'e'd by the body’s rotation about its 
axis). The C 22 te'rm elescribes the ellipsoielal nature 
(the ecpiatorial ellipticity) of a body, and the J;j co- 
efficient models the nonspherical mass distribution 
between the northern and southern hemispheres (the 
body’s “pear-like” shape). Higher degree gravita- 
tional terms describe more localized distributions of 
nuiss. Figure 4 (from ref. 22) illustrates the equi- 
potential surfaces for the low-degree zonal harmon- 


Figure 4. Geometrical sliape of Legendre polynomials cor- 
responding to ecpiipotential surfaces for zonal liarmoii- 
ics. The surfaces shown hen' an' for positiv(' values of 
J coefficients. (From ref. 22.) 

ics. The degree of the term determines tlu' number of 
lobes of its equii)ot(uitial surface. As previously men- 
tioned, the zonal cquipotential surfaces are indepen- 
dent of longitude and model axially symmetric po- 
tentials. The cquipotential surfaces associated with 
the sectorial and tesseral harmonics (figs. 5 and 6) 
are functions of both latitude and longitude (the 
z-axis shown in the figures corresponds to the spin 
axis of tlu^ Moon, and the positiv<' :/:-axis is aligned 
with the Moon- Earth direction). For sectorial and 
tesseral harmonics, the order of the term represents 
the number of lobes of a horizontal cross section 
of its cquipotential surface (corresj)onding to cos 7;/ A 
and sinmA). 

Values of the gravitational coefficients are deter- 
mined from radar tracking data of lunar satellites 
and laser-ranging data measurements. (S('c rc'f. 23.) 
With the Doppler method, the radial velocity of the 
satellite (the velocity component along the line of 
sight from the tracking antenna to the spacecraft) is 
determined from the difference in frequency l)(^twecn 
the signal emitted by the satellite and the signal re- 
ceived by Earth tracking stations. From tlu'se ve- 
locity measurements, the radial accek'ration can be 
determined. By attributing this acceleration of the 
satellite to the gravitational field of the Moon, the 
potential field (and therefore the value of the har- 
monic coefficients) can be derived. These coefficients 
arc independent of the orbit of the spacecraft lur- 
ing tracked, because they describe the fluctuations 
in the gravitational field. However, in practice, the 
estimated coefficients are dependent on the altitudes 
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and inclinations of the satellite data used to derive 
the model. This dependence is a result of the limited 
geographic distribution of the measurement set used 
to compute the coefficients. All effects of the lunar 
gravity field have not been included in a coefficient 
model; a set of satellite data covering all regions of 
the Moon is as yet unavailable. Also, the number 
of terms included in the gravitational model (the 
maximum value of the index ji) also influences the 
estimated values of the coefficients. 

Lunar laser-ranging data arc obtained by mea- 
suring the round-trip time of a laser pulse between 
an Earth observatory and a retroreflector array on 
the surface of the Moon. These measurements allow 
the Earth-Moon distance to be determined within 
fractions of a meter, precise enough to monitor the 
physical librations of the Moon. By observing varia- 
tions in the lunar rotation, values for the low-degree 
harmonics can be determined. (See ref. 24.) 

Limitations of existing models. There are 
three major obstacles to making an accurate de- 
termination of the Moon’s gravitational field. (See 
ref. 25.) First, ^is previously mentioned, the various 
surface and interior features make the Moon a very 
coinplicated gravitational object, and therefore diffi- 
cult to model with a simple mathematical represen- 
tation. Further, gravitational models of high order 
and degree may not be appropriate for orbital life- 
time predictions at the present time, because there 
is still uncertainty in the values of even some of the 
lower coefficients, such as J;^. These uncertainties 
may overshadow any attempt to adequately model 
localized variations in the Moon’s gravitational field 
with higher coefficient gravity models. 

Second, the range of inclinations and other orbital 
parameters for which tracking data exist is quite lim- 
ited. Ample tracking data exist for near-equatorial 
orbits. However, for mi<llatitude and near-polar in- 
clinations, little tracking data are available. Since 
efforts have been made to include only orbital data 
with trajectories free of propulsive maneuvers (influ- 
enced only by the force of gravity), the amount of 
tracking data available is further limited (ref. 26). 
High-altitudc orbits are ideal for determination of 
the low-degree harmonic coefficients, because these 
orbits are not sensitive to the effects of high-degree 
harmonics. (See ref. 27.) Howwer, most of the lunar 
satellite tracking data available are for low- altitude 
orbits. The last U.S. satellites to orbit the Moon 
were the Explorer probes in the mid-1970’s. Little 
new information about the Moon’s gravitational field 
has been obtained since then. 


The inability of Earth-based stations to obtain 
satellite tracking data as the spacecraft passes he- 
hind the Moon is the third obstacle to developing 
an accurate lunar gravitational model. Because of 
libration of the Moon’s orbit, 41 percent of the lunar 
surface is never visible from the Earth (ref. 28). As 
a result of this limitation in tracking data, gravita- 
tional effects on the far side of the Moon cannot be 
practically determined. The absence of a description 
of the gravitational field on the far side of the Moon 
remains as the greatest limitation to making accurate 
predictions of orbital lifetimes for lunar satellites. 

Current models available. Several lunar grav- 
ity models have been constructed, based primar- 
ily on satellite tracking data, with the additional 
use of lunar laser-ranging observations and mass- 
concentration models. Differences in the values of 
the gravitational coefficients of these models occur 
because of the selection of different satellite trac:k- 
ing data and the different methods used to process 
these data. (See tables 1 and 2 .) Coverage from 
the Apollo satellites is limited to about 20 percent 
of the lunar surface. (See ref. 29.) Tracking data 
from Lunar Orbitcr 4 and 5 are commonly used to 
provide high-inclination information, while Ai)ollo 15 
and 16 subsatellitc tracking data an' used to {)rovide 
low-inclination information. 

Six gravitational models wert' recently compan'd 
by the Jet Propulsion Laboratory for the Lunar Ob- 
server mission study. (See ref. 30.) The Liu-Laing 
model (ref. 31), constructed primarily from Doppler 
data of the five Lunar Orbiter satellites, is an 8 x 8 
(n X m denotes nth degree, mth order, see eq. (3)) 
model with additional zonal harmonics up to de- 
gree 15. This model differs from other models in 
that the values of the high-degree zonal harmonics 
are larger than the value of the low-degree zonal 
harmonics. The Ferrari 5x5 model (ref. 23), d('- 
veloped from 9 days of Lunar Orbiter 4 data and 
approximately 2200 lunar laser- ranging obs('r vat ions, 
is an attempt to accurately determine the values 
of the low-degree coefficients. A Ferrari 16 x 16 
model (ref. 26) was also developed with Doppler 
measurements from Lunar Orbiter 5 and Apollo 15 
and 16 subsatellites. The Bills-Ferrari 16 x 16 model 
(ref. 32) used the data sets of the Ferrari 5 x 5 
and 16 x 16 models along with a model of approx- 
imately 600 rnascons. Akim (ref. 33) developed a 
4x4 model with zonal coefficients J 5 , J{j, and Jj, 
based on the Soviet Luna spacecraft. A Sagitov 
16 x 16 model (ref. 34) also exists, based on data 
used for the Ferrari 5x5 model, the Bills-Ferrari 
model, and the Akim model, with additional data 
from Apollo and a rnascon model. As mentioned 
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Bills-Ferrari 16 x 16 

Sagitov 16x16 

Akim-Vlasova 16 x 16 

Ferrari 16x16 

Ferrari 5 x 5 

Liu-Laing 15x8 



Figure 7. Comparison of six gravitational fields used in .IPL 

Lunar Observer study. Initial 100-kin circular ]X)lar orbit. 

(From ref. 3b.) 

previously, all six of these models are limitetl by a 
lack of tracking data for spacecraft on the far side of 
the Moon and by a scarcity of low-altitude, high- 
inclination data. Although these models generate 
similar results for equatorial orbits, the results pre- 
dicted by these models vary significantly for an ini- 
tial 100-km-altitude circular polar orbit. (See fig. 7, 
taken from ref. 35.) 

Much debate currently remains as to which grav- 
itational model provides the best predictions of or- 
bital lifetimes for lunar satellites. In the present 
analysis, the consequences of choosing a particular 
gravitational model are not emphasized; the focus is 
centered on the selection of the initial orbital param- 
eters. The Ferrari 5x5 model (ref. 23) was chosen 
for this study for several reasons. First, the computer 
program used in this analysis (appendix B) is limited 
to gravity models of degree and order no larger than 
eight. Second, this model provided a '^vorst case” 
scenario for polar orbits, because this gravitational 
model generated the shortest orbital lifetime predic- 
tions in the JPL study. Therefore, this model may 
overestimate the rate of decay in perilune altitude. 
However, for mission design purposes, overestimation 
of the effects due to the MoorFs gravitational field is 


more desirable than underestimation. A low-degree 
model was also desired to avoid addressing the errors 
and effects of the higher degree harmonics. Although 
the low-degree model is incapable of modeling local- 
ized gravitational anomalies, it provides a global de- 
scription of the lunar gravity field that is adequate 
for performing lifetime studies. 

Analysis 

The goal in this analysis is to illustrate the im- 
pact of the Moon’s gravitational perturbations on 
the lifetimes of low-altitude, near-circular parking or- 
bits, and to identify orbits favorable for lunar out- 
post missions (orbits that have the longest orbital 
lifetimes). A method is provided for mission plan- 
ners to incorporate nonspherical gravitational effects 
in preliminary lunar analysis studies. 

Effects of External Forces 

The Lagrange planetary equations presented in 
equations (4a) to (4f) describe the changes in a satel- 
lite orbit due to forces that may be expressed in terms 
of a pot(mtial. The equations in appendix A have 
assumed a potential solely as a result of the lunar 
gravitational field. To give a complete description 
of the satellite motion, all forces acting on the satel- 
lite must be taken into account. The absence of a 
significant atmosphere on the Moon eliminates any 
need to consider drag effects. Forces that could affect 
the satellite motion include solar radiation pressure 
and perturbations from other bodies, particularly the 
Earth and the Sun. 

The influence of solar radiation pressure on a low- 
altitude lunar satellite is examined in this report. 
(See appendix B.) The resulting acceleration of a 
satellite is a function of its mass, cross-sectional area, 
and reflec'tivity coefficient. A mass of 10.2 metric 
tons and a cross-sectional area of 75 nr were chosen 
as values for a lunar excursion vehicle; also, the 
vehicle' was assumed to be a perfect reflector. Solar 
radiation pressure has the greatest effect on low- 
inclination orbits, although these effects were small. 
For a near- circular orbit with an initial altitude 
of 300 km, solar radiation causexl a change in the 
perilune altitude of less than 1 km over a period 
of 180 days. Unless a satellite has a large cross- 
sectional area and small mass, solar radiation tffiects 
are negligible. 

Although both the Earth and the Sun contribute 
perturl)ation effects, the effects due to the Earth arc 
about 170 times larger than the effects due to the 
Sun. (See ref. 36.) The average effects of a third body 
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oil the eccentricity (the main parameter in orbital 
lifetime studies) can be expressed as 
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where the bars denote variable's that have been aver- 
aged over the period of the third body about the 
central body and over the period of the satellite. 
(See ref. 37.) Tlu' variable u/j represents the mean 
angular motion of the third body about the central 
body. From equation (5), it is apparent that third- 
body perturbations have the greatest effect on eccen- 
tricity for mid-eccentricity (maximized for c — 0.707) 
and high-inclination satellite orbits. The average 
value of Lj also strongly influences the effects of the 
third-body perturbations. 

Using the program LUNLIFT] with the Ferrari 
5x5 gravity model, a near-circular, initial 300-km, 
perilune-altitude polar orbit had a 20-kin-higher per- 
ilune altitude after 180 days when third-body per- 
turbation effects were included. (See fig. 8.) For 
an initial 100-km perilune-altitude polar orbit, the 
orbital lifetime Wcus 157 days when third- body ef- 
f('cts were included; the lifetime was only 144 days 
wlu'ii Earth-Sun {)erturbation effects were lu'glected 
(fig. 9). Although the effects of third-body pertur- 
bations on orbital lifetinu's may not be negligible for 
all orbits, these effects are small for most of the or- 
bits addressed in this investigation (bec'ause of the re- 
st rid ion of near-circular parking orbits). Therefore, 
third-body pc'rturbation effects were not included in 
this analysis. 


Earth-Sun perturbations included 



Time, days 


Figure’ 8, Farth-Sun perturbation effects for 300-krii initial 
pcTiluno altitu(k‘ polar orbit (using Ferrari model). (Initial 
conditions: i = 90”. il = 0'\ uj ~ 225°.) 


Earth-vSun perturbations included 



Figun; 9. F^artli-Sim iH'rturbatum effects for 100-kin initial 
jx'riluiK’ altitude polar orbit (using F(‘rrari model). (Initial 
conditions; i — 90°, il -- 0°, u.’ - 225 ’.) 

Use of Orbital Elements for Lifetime 
Studies 

The Lagrange planetary equations of iiKffdon 
(eqs. (4a) to (4f)) are used for orbital lifetime stud- 
ies, because they enable integration of the classical 
orbital elements (u, c, /, fh cj, and A/, see figs. 10 
and II). The mean anomaly A/, dcdiiu'd as th(' rru'an 
angular motion multiplied by th(' tintt’ since peri- 
center passage, is related to the true anomaly u. The 
advantage of utilizing orbital elements as opposed to 
Cartesian coordinates (dchning the radius and v(v 
locity vectors of the satellite as functions of y, 
and z) is that v(Ty large' integration time steps can 
be used, because^ the orbital elements (except for A/) 
change vx'ry little over coiise'ciitivc' orbits compared 
with a Cartesian de'scuiptioii. Integration over large' 
time' ste'ps require's that the' short -perie^ei e'ffbcts be' 
averageel fe)r t he disturbing function in the' Lagrange' 
c^quations. Time steps e)f 1 day we're use'd to geiK'rate 
the results in this analysis. To verify that 1-day time' 
steeps were not too largo to give' inae-euiratc re.'sults, 
test crises were run with the program LUNLIFF] with 
10-sec time stops. The re^sults generateel lyy the'st' 
two different step sizes we're' nearly ieleritie’al. within 
0.01 perce'iit. 

The disadvantage e)f using a classical orbital- 
eleuncnt elescriptie)ii of the' e'epiatiems of motie)ii in- 
volves the inability te) directly addre'ss certain type's 
of orbits (e.g.. circular, eeiuatorial), be'causc some 
of the elements may not be elefine'd. The' classi- 
cal orbital-clement description eloe's not pose? a sig- 
nificant restriction in this stuely. as near-circular 
(e = 0.05) and near-equatorial (/ ~ 1°, aiiel 179°) 
orbits were addressed. (There is a se»t of orbital 


10 




Fifz;urt‘ 11. Orbital geoinotry. 

elciiients, referred to as equinoctial elernents, that 
can address circular and equatorial orbits. (See 
refs. 30 and 38.) However, (equinoctial elements arc 
utilized much less frequently than classical orbital 
elements.) 

Only three of the six orbital elements were moni- 
tored in this analysis. Since the mean anomaly M is 
used to describe the position of the satellite within 
the orbit, and not the shape or position of the orbit 
itself, it need not be considered in determining orbital 
lifetiiiK's. The semimajor axis experiences only short- 
period variations, so as it is integrated ovor the mean 
anomaly for one period, its average rat(' of change is 
zero. Short- period effects in the orbital elements are 
observed with the program ASAP (Artificial Satellite 
Analysis Program, ref. 39). The amplitude of short- 
term variations is small (on the order of 1 km (figs. 12 
and 13) and has little effect on orbital lifetime pre- 
dictions. Hence, the semimajor axis is treated as 


Short-period effects included 

(program ASAP) 


Short-period effects not included 

(oroeram LUNLIFE) 



Time, days 

Figure 12. Slioit-jH'riod ilFivts in st'iniiimjur axis. (Initial 
conditions: i = 90b H ~ 0”, a.' = 0°. fip — 300 kin: gravity 
inod(9: 8x8 truncation of Hills-F('rrari niod('L) 

Short-period effects included 

(program ASAP) 

Short-period effects not included 



Time, days 

Figure 13. Short-ix'riod effects in pcudlune altitude'. (Initial 
conditions: i = 90°, 1] = 0°, u; = 0°, ftp - 300 km; gravity 
model: 8x8 truncation of Bills- Ferrari model.) 


a constant in this analysis. Inclination is cyclical 
and has variations typically on the order of 1° over 
a 180-day period. (See fig. 14.) As with the stmii- 
rnajor axis, the inclination can be trc'ated as a con- 
stant with little loss of accuracy (although chang(\s 
in inclination were included in this analysis). The 
remaining orbital elements that change significantly 
over the time periods studied in this analysis (other 
than the mean anomaly) are eccentricity, longitude 
of ascending node, and argument of perilune. 
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Kigufo 11. Typical inclination versus time {)lot. (Initial 

conditions: / ~ 1°, — 0^^, u; = 0^, hp = 100 krn.) 

The key parameter in orbital lifetime studies is 
the perilune altitude (or equivalently, the eccentric- 
ity, sinc'c the perilune altitude is a function of the 
ecc(mtricity). The life of the satellite orbit is termi- 
nated when the magnitude of the perilune radius of 
the satellite falls below the value of th(' average radius 
of the Moon; the satellite then impacts the surface. 
Th(' perilune radius is given as: 

'> = o(l - (’) (6) 

Simple differentiation, iissurning the stMiiimajor axis a 
remains c’oiistant, yields 

A/‘p = —a Af’ (7) 

The rate of change in peidlune altitude is directly 
proportional to the rate of change in eccentricity. 
Near-z('ro or i negative rates of change (orbit becomes 
mor(' cii’cular) in ec'ctmtricity enable orbital lifetimes 
to be maximiz(xl- Since the rate of change in ec- 
centricity is a function of initial inclination, longi- 
tude of as(‘emling node, and argument of perilune 
(eqs. (A2b), (A3b), and (A5b)), these initial orbital 
paranud(n‘s directly affect tlu' orbital lihdime. Per- 
ilune altitude is also indirectly affected by the rates of 
chang(' in longitude of ascending node and argument 
of perilune. 

The initial orbital parameters were chosen to in- 
clude many types of parking orbits that would be 
of interest to mission planners. The initial eccen- 
tricity was fixed at 0.05, small enough for the orbit 
to be considered nearly circular, yet large enough to 
avoid any singularities in the computer program asso- 
ciated with circular (zero-eccentricity) orbits. Circu- 
lar parking orbits are attractive to mission planners, 


because they simplify phasing requirements for de- 
scent to and ascent from th(' lunar surface, and sim- 
plify phasing re(iuirements for an Earth-return burn 
(although higher AP values are required to achieve 
circular parking orbits). Circular parking orbits are 
also generally preferred for mapping missions of the 
Moon's surface. For this study, initial perilune alti- 
tudes of 100 km and 300 km, referenced to the mean 
value of the radius of the Moon, were selected. The 
specific inclination value of a parking orbit selected 
for a particular mission will be heavily influenced by 
the desired mission objectives and the location of lu- 
nar landing sites. In this analysis, both direct and 
retrograde orbits are addressed. Lifetime results were 
first generated by varying initial values of longitude 
of ascending node and argument of perilune over 360° 
for a specific initial value of inclination. Since orbital 
lifetimes were weakly influenced by initial values of 
longitude of ascending node (small effect of the sec- 
torial and tesseral harmonics on lifetimes), lifetime 
results were generated by varying initial values of in- 
clination and argument of perilune for fixed vahu's of 
initial longitude of ascending node. 

Results and Discussion 

Development of a Simplified Gravitational 

Model 

The need to analyze orbital lifetimes for a large 
number of initial orbital parameters motivat(xl the 
formulation of a simplified gravitational field model. 
Although the program LUNLIFE was capable of gen- 
erating the lifetiiiK' of a given orbit relatively cpiickly 
with the Ferrari 5x5 gravitational model, excessive 
computer time (on the order of 100 hours) would 
have been required to geru'rate the large number of 
results desired in this analysis. For example, to gen- 
erate contour plots of lifetimes versus initial argu- 
ment of perilune and initial longitudes of asce'iiding 
node to a resolution of a 5° by 5° grid, the eval- 
uation of approximately 5000 orbital lifetimes was 
necessary. By neglecting many of the coefficients in 
the Ferrari 5x5 model and using a simplified in- 
tegration technique (a single-step integration of th(^ 
hrst-ordcr equations of the orbital rates instead of a 
multiple-step numerical integration of the complete 
equations as given in ref. 40), a mor(‘ (efficient method 
was d(wis(xi to predict orbital lifetimes. This method 
reduced the required computation time by approxi- 
mately two orders of magnitude w'ithout apprecia- 
bly sacrificing accuracy (less than 5 percent) in tlu' 
prediction of orbital lifetimes. 

Previous studies have attempted to approximate 
the MooiTs gravitational field by using only a few 
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harmonic coefficients. A four-coefficient lunar grav- 
ity model was developed in the 1960’s (based on Lu- 
nar Orbiter tracking data) for Apollo mission control. 
(Sec ref. 3.) The formulation of this simplified model 
was necessitated by the restrictions of onboard com- 
puting capabilities at that time. The coefficients J 2 , 

^ 22 - C 31 were selected for this model because 

it was determined that they had the greatest effect 
on the orbital elements. The model was limited in its 
design to accommodate only initial orbital elements 
similar to the Lunar Orbiter missions. A similar sim- 
plified model that was developed for this analysis has 
been validated for the entire range of values of initial 
inclination, argument of perilune, and longitude of 
ascending node. 

A reduction in the number of computations also 
allowed for the possibility of implementing optimiza- 
tion routines to maximize orbital lifetimes. An 
attempt was initially made to utilize optimization 
methods for selecting the initial orbital elements that 
yielded maximum lifetimes. However, this approach 
was deemed infeasible because of the highly non- 
linear nature of the problem. The optimization re- 
sults were also dependent on the initial guess, be- 
cause the orbital lifetimes objective function w^is 
characterized by many local maximum and minimum 
values. 

Several issues were taken into account when se- 
lecting the coefficients for the simplified lunar gravi- 
tational model. First, an effort was made to retain all 
the terms used in the Apollo mission control model. 
Second, the numerical values of the coefficients were 
examined to decide which ones could be eliminated. 
(See table 1 .) An attempt was made to kcx^p as many 
zonal terms as possible, because they are responsi- 
ble for secular effects and arc therefore significant 
in orbital lifetime studies (the J 4 term was eventu- 
ally discarded because of its small magnitude). Fi- 
nally, the lower order coefficients were assumed (and 
later verified) to have a greater effect on the variation 
in the orbital elements than the higher order terms. 
A sample case was selected and analyzed with the 
full Ferrari 5x5 model. Individual coefficients were 
discarded, and the sample case was evaluated and 
compared with the full model after each term was 
eliminated. All the Snjn tc'rms were quickly elimi- 
nated because of their small value. The higher order 
Cjim teams were eliminated for the same reason. The 
final coefficients selected for the simplified model con- 
sisted of the same terms as the Apollo mission control 
model (J 2 , J 3 , C‘ 22 ^ and C 31 ), with the addition of 
the J 5 term. The values of these coefficients were 
adopted directly from the Ferrari 5x5 model. An 
attempt to eliminate the J 5 term drastically changed 


the lifetime results. This change indicated the need 
to include the J 5 term when performing orbital life- 
time predictions for long lunar stay times with the 
Ferrari 5 x 5 model. 

The orbital lifetime results generated with the 
simplified five- coefficient model corresponded closely 
to the results generated with the Ferrari 5 x 5 model 
(appendix C). Comparisons of the plots of the or- 
bital elements versus time for the two models illus- 
trate the higher resolution prediction capability of 
Ferrari's model. (See figs. 15 and 16.) The simpli- 
fied model contains fewer coefficients to contribute 
periodic effects and thus simplifies the shapes of the 
graphs. 



P^igun' 15. Coniparison of ]H'riliine altitudt' for simplified and 
Ferrari models. (Initial conditions: i = 120°. H = 0°, 

u; = 225°, h,> = 100 km.) 



F^igure 16. Comparison of argument of perilune altitude 
for simplified and Ferrari models. (Initial conditions; 
i ^ 120°, Q = 0°, a; = 225°, hp = 100 km.) 

The simplified five-coefficient gravity model was 
an important tool for generating results quickly and 
for simplifying the interpretation of the final results. 
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Since first-order analytical equations for the orbital 
rates are readily available in the literature (eqs. (Al) 
to (A5)), a simple single-step integration scheme can 
be used with these equations. For example, the 
iiiti'gration of eccentricity is pi'rfornu'd by using the 
('quation 

~ CM + At (8) 

where represents the finite-integration time step 
and ^ represents the superposition of the eccentric- 
ity rates of change contributed by each gravitational 
coefficient (assumed to be constant over the time 
step At), Generating the major effects of tlu‘ Ferrari 
5 X 5 model with only a few coefficients also aided 
in the interpretation of the results. Since first-order 
ecluations existed for each coefficient in the simpli- 
fied model, analytical methods could Ix' applied to 
(’xplaiii the ('ffects of individual coefficicuits and the 
sujx'rposit ion of tlu'se effects. 

Effects of Individual Gravitational 

Coefficients 

d'lu^ gravitational coefficients in the simplified 
gravitational model giv(' ris(^ to both s('cula,r and [)e- 
riodic effects. To first-t)nh'r analysis, sc'cular effects 
arc’ prodiu'ed only by the zonal harmonic tca'iiis. Both 
t,h(‘ zonal and Jr^) and the’ off-diagonal terms 

(Cb'i- ^ai) ill fids model are responsible’ for pc’riodic 
(dfc'C'ts. Cxxdficients in the lunar gravitational model 
give rise to short-period, medium- period, and long- 
[)('riod variations. The short-pc’riod (’ffc'cts (on the 
ordc’r of the’ satc’llite's orbital period) wc’rc’ not mod- 
(’led in this analysis. Tlmse efh’cts can bc' nc’gk’ctc’d 
to first approximation, because tlx’ rate’s of change 
of thc'sc' c’ffects cjften average’ te) zerey as is the Ciuse 
for the’ semimajor axis. The meMlium-pcriexl terms, 
e’ommemly refe’rnxl to as the lunar /n-daily terms, 
are’ on the ordc'r of a lunar day, e)r 27. d Earth days. 
The’ infliu’ue-e of this pc’riod can be easily seen in 
the' plots of the orbital elements versus time (figs. 17 
and 18). A se'ctorial or tc'sse’ral harme^nic term of 
oreler m produce's nu'dium-pe'riod e’ffects with a pe- 
riod of about 27.3/m days (the exact period varie’s 
slightly as a re'sult of pre'ce'ssion of the’ longitude of 
ascending node). Higher order harmonic tc’rms have 
k'ss iiiHueiu'e on long-term orbital lifetime studies, 
lx’e*aus(' short -period effects are integrate’d out more 
(luieTly than long-period effects. Long-period effects, 
on tlu’ order of a year, can also be obse’rved in the 
variations of some’ of the orbital e’leme’iits. Thc’se ef- 
fects, contributc'd by the zonal harmonics, arise from 
the vsec'ular variation of u;. The long-period effects 
are important for lifetime studies, be’cause they pro- 
duce the largest contribution in pe’rilune’-altitude de- 



Figure 17. M(‘dium-p(‘rio(l (’fleets in argmiK’nl of jx’rilurH’ 
for F(’rrari model. (Initi^d e'ondit ions; i = 179'', il - 0 \ 
uu' = hff — 100 km.) 



Figure IS. M(’dium-])(‘i iod effeets in perilum' for teirari 

model. (hiitial conditions: i - 179'\ H ~ (F. — (F. 

fij) = 100 km.) 

cay rates over long pe'ritxls of time and behave' as 
secular effects for orbits with lifetime?s of k’ss t han a 
year. 

The rate’s of ediangc in c. SL an el uj are eiire’edly 
pro])ortional to the value's of the gravitational co- 
efficients, although not every coefficient causes a 
change in every e)rbital ekunent. (Se'c table 3.) 
The J 2 te’rni is the elominant term for de’te’rmining 
the rate e>f change of longituele e>f asceneling node. 
This term and the C 31 term also infiue'iice’ the^ rate 
of change in argument e>f perihme. The' 
and C 31 terms clireedly affe’ct the rate of change’ in 
eccc'iitricity anel, therefore, perihme altitude’. The 
final term in the simplifienl model, 622 - primarily 
influences the orbital inclination. 
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Orbital Lifetime Study 

A single-step (Euler) integration technique^ was 
used with the simplified fivc-coelficient gravity model 
to generate all results in this analysis. Onee initial 
orbital elements were specified, a particular orbit was 
analyzed for a period of 365 days, unless the lifetiiiK' 
of the orbit expired in less than a year (in which 
case the analysis was halted immediately). If the 
lifetime was less than 365 days, the actual value of 
the lifetime is pres('iited. If the lifc'time exceeded 
365 days, the minimum i)erilune altitude attained 
during the 365-day period is j)r('s(^nted. The program 
LUNLIFE, which uses tlu^ complete Ferrari 5 x 5 
gravity model, was periodi(‘ally eiiii)loyed to provide 
validation for lifetime predictions of individual orbits. 

The terms J 3 . J 5 , and were the only co- 
efficients considered that dir{H:tly aifi'ct rates of ec- 
centricity and, therefore, orl)itaI lifetimes. Unlike 
the zonal harmonic coefficients, the eccentricity rate 
for is a function of initial longitude of ascend- 
ing node and introduces medium-pc'riod effects with 
a period of about 27.3 days. This tc^rm has a dra- 
matic effect on rates of eccentricity (but. tlu'se are 
only localizcnl effects), and giv(^s rise' to step- function 
behavior in plots of lifetinu^ (fig. 19). This figure* di- 
rectly illustrates the effect that C:^\ has on orbital 
lifetimes. 



Figure 19. ion Ijclmvior in lik'tinu's j>rofil(’ (hie to 

Un ttTii) (/ = 55""). 

The st(q)-function sliaiH* of the* lifetiim* plot is 
a result of the* superposition in eccentricity rate\s 
between anel the zonal coeffide'nts. The* C:u ten'iii 
e'ontribute'S large m(*dium-pe*riod e*ffe*('ts, but the ,/;$ 
anel J 5 terms cenitrilmte long-i)e*rioel effects. When 
the decay rate is of the same* sign as the* ze)nal 
ele*cay rate, the insult is a sharp drop or rise in e)rbital 



(a) (\>ntour ])lot. 



(b) T]ir('('-{liin(‘usioiial plot. 

P'igurt* 20. Minimum ix’rilune all iiink* {l(’])(‘ii(l(nK“(' on initial U 
and u.' for low inclinations (/ = 3"). lOO-km initial pjn'iluiK’ 
altitiuk’ ov(‘r a 3f)5-day period. 

lifetimes. Whe*n the ("yq de*cay rate is of the oppeisite* 
sign e)f the* ze)]ial de*cay rate, the total ele*cay rate* is 
nc*arly ze*rex tins low rate insults in plat(*aus in the* 
pk)t e)f or!)ital lifetime's. 

The eff'eed e)f initial longituele e)f asce*neling iioele 
on orbital ele*e‘a.y rate's was significant. Althe)ugh 
orbital Ihetime insults ajtjie'a.r te) be (le*i)e*nde*nt em 
the* initial longitude* e>f asenneling noele for ne*ar- 
equatorial orbits (le*ss than 1 (U inclination), the*y 
are actually depenelent em the* sum of longitude of 
asceneling noele anel argument e)f pe*rilune. Figure* 20 
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shows the variation in minimum porilune altitude 
as a function of initial argument of perilune and 
longitude of asctuiding node for a 100 -km initial 
perilune- altitude orbit with a 3° inclination. As 
shown in this figure, various initial orbital conditions 
are not described by independent values of and cj, 
but rather by values of -h u;. For this reason, 
orbital lifetimes arc nearly the same for orbits along 
lines of constant + (See fig. 20(a).) Lifetime 
results were subsecpiently generated by varying the 
initial argument of perilune and inclination over the 
entire range of possible values. These results are 
independent of the initial longitude of ascending 
node, although slight differences in the results will 
arise, depending on the value of 12 . 

Figures 21(a) to 21(c) show the instaiitaneous 
pi'rilune altitude decay rates for 17 = 0°, 45°, and 90°. 
Although these plots exhibit different characteristics 
for different values of 12 . similar orbital lifetime re- 
sults will occur regardless of the value of 12 , because 
the effects due to 17 will average out to zero after 
a period of 27.3 days. The plots only illustrate the 
instantaneous decay rates. To generate valid orbital 
lif(4iiii(^ predictions, these decay rates must be in- 
tegrated o\x'r time. Although the inclination varies 
little, the value of argument of perilune may change 
greatly. An orbit with a given set of initial condi- 
tions such that the altitude decay rate is large may 
evolve over time to a set of conditions for which the 
decay rate is extremely small (primarily caused by 
a variation in argument of perilune); therefore, inte- 
grated lifetimes longer than those inferred from the 
figures may be exhibited. Nevertheless, the figures of 
perilune altituch' decay rates are useful in illustrat- 
ing the magnitude of decay that a lunar satellite can 
experience in a day. 

An integration of the orbital element rates for or- 
l)its of various initial conditions was utilized to pro- 
duce figure 22. The large white regions occur where 
the initial orbital parameters had a lifetime in excess 
of 3()5 days. Tlu' r(\sults display a symmetry about 
i =z 90°^ which is violated slightly by the presence of 
a cosine i term in some of the orbital element rate 
ecpiations. The closely spaced contour lines in the 
figure are the result of localized effects due to the 
C;ji coefficient. 

Fiv(^ bands of sliort-lifetime orbits occur for spe- 
cific values of inclination. The appearance of these 
bands can be explained by examining the perilune al- 
titude decay rates (directly proportional to the nega- 
tivi' of the eccentricity rates, see eq. (7)) for J 3 and J 5 
cus a function of inclination. (See fig. 23.) Figure 23 





Figure 21. Perilune-altiludc' decay rates in kiii/day vitsus i 
and LO. 
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Initial inclination, i, deg 


suit of contributions from the J3 term. Also, bands 
of long- lifetime orbits correspond to zero values of 
the eccentricity rates for the J5 (coefficient [i = 0°, 
40.09^ 73.43°, 106.57°, 139.91°, and 180°). 

The J5 coefficient is the main driver in the 
Ferrari model for predicting lifetime effects. Further 
evidence of the Jr^-term dominance is shown in fig- 
ure 24. This plot of orbital lifetimes cis a function 
of initial argument of perilune and initial inclination 
was generated with the J3 and C;p coefficients set 
ecpial to zero. This figure is similar to figure 22, which 
included effects from J3 and C31, and shows that the 
,/5 term contributes the dominant effect to perilune 
altitude decay, although the effects from J3 and C31 
arc certainly not negligible. (See figs. 19 and 23.) 


Figure 22. Contour plot of orbital lifetimes in days for initial 
longitude of fUscending node of 0°. Larg(" white regions 
represent lifetimes longta* tiian 1 year. 


Superposition of J3 and rates of change 

in perilune altitude 


Rate of change in perilune altitude due to J3 

Rate of change in perilune altitude due to J5 



Figure 23. Dependence' of rate of change in perilune altitude 
on inclination for J:\ and Jr, terms. 100-krn initial perilune 
altitude'. 

was generated by multiplying equations (A3b) 
and (A2b) by —a and plotting them, with the initial 
orbital conditions and a; — 0. (All terms proportional 
to in the ^ equation for J5 were neglected, bev 
cause their contribution is small for near-circular or- 
bits.) The difference in the amplitudes of the wave- 
forms illustrates that contributions to the decay rate 
are much larger for J5 than for J3 in the Ferrari 
model. The bands of short-lifetime plots roughly 
correspond to the minimum and maximum peaks 
of the J5 eccentricity rate {i = 19.42°, 56.14°, 90°, 
123.86°, and 160.58°) and are slightly offset as a re- 



30 45 90 120 150 180 

Initial inclination, i, deg 


Figiin' 24. Contour plot of lifetimes iii days for initial lon- 
gitude of ascending node of O'C lOO-km initial periluiu' 
altitude; effects of J\\ and are neglected. 


Results were also generated for initial perilune 
altitudes of 300 km. Figure 25 is a plot of minimum 
perilune altitude as a function of initial argument of 
perilune and inclination. The regions at inclination 
values of approximately 55°, 90°, and 125° represent 
areas where the orbit crashes in less than a year. 
Although the orbital decay for the higher altitude 
orbit is less than for the 100- km case, the pattern 
for the 300-km case is very similar. This analysis 
focuses primarily on the 100-km initial altitude case, 
with the assumption that results for this case can be 
correlated with the 300-km initial altitude case. 

Sensitivity Studies — An Assessment of 

Uncertainty in Gravitational Field 

As previously mentioned, many problems have 
been encountered in developing an accurate descrip- 
tion of the Moon’s gravitational field. Figure 7 illus- 
trates the discrepancies among the current available 
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Figun' 2"). C’uutour plot of iiiiiiiiiiuui piTiluiu' aUitudf in kin 

()V(M‘ a dlin-day poriod. dOO-kin initial j)oriluiR‘ aJlilude. 

iiiodids. Two (piestions naturally arise*: How aceai- 
ratedy are the eoeffieients in the F"(Trari 5x5 model 
known (how much error is assotdated with each co- 
(dfiedimt?) and is the Ferrari model (or any other 
mod(d) an accurate description of tlie actual lunar 
gravitational held? The second question is beyond 
tlu* sc'ojK* of this investigation, but tin* first cpiestion 
(‘an Ik* addre'ssed by examining the standard devi- 
ations in the values of the gravitational coeffici(*nts 
statc'd for the Ferrari 5 x 5 model. Mor(' importantly, 
tlu* possible (uror in orbital lifetime preditdions as 
a rt'sult of th(‘ imcindainty in tlie gravitational co- 
(dhtdemts can be addrt'ssed with the* use* of sensitivity 
coedhedemts. 

The* tormulation of the sim])lifie*d gravitational 
moded illustrate'd that five* coe'fhedents we*re adeejuate* 
to (*xj)lain the effects due to the Ferrari model on 
an orbiting lunar satellite. H(*nce, a good (*stiinate 
of the* ('irors associatexl with the unce*rtainty in the 
value* of tlie* coefheients for the Ferrari model can 
be* e)btaine*d by (xmside*ring only the uncertainti(*s of 
the hve coefficieeits inedude'd in the simpliheMl model; 
these five e'oe'fhciemts have* the* gre*ate‘st influence on 
the* orbital e*k*me*nt rate's. Fdirthermore, since only 
thre*e of tlie'se (coefficients (Jp J 5 , and (7;ji) dire'ctly 
affect the* ('(‘('('ntricity i-ate*s (the J? arid Cr) terms 
indire*e‘tly affect the* e'(*(‘e*ntri(dty rate's by causing 
change's in the otlu'r orbital ele*nients), the e*rror 
in orbital lifetime' pre'dictions associated with the 
unc(*rtainty of the* coefficients in the Ferrari moded 
can be approximate'd by addre'ssing the uncertainty 
sede'ly in the*se* thre'e harmonic terms. 

The* sensitivity coefficients for each gravity co- 
effie'ient are shown in figures 2C to 28. The sensitivity 



Figure' 2G. Contour plot of s(*iisitivity for J:\ in units ed' 
kin/day x 10‘^ = 0°) for ftp = 100 kin. 



Figure* 27. Contour plot of se'iisitivity for Jr> in units of 
kiu/(iay X 10^ (U — O'd for ftp = 100 km. 

coe'ffiedemts we*re* obtained by differentiating the* rate* 
e)f change in perilune altitude (the first-order ana- 
lytical equation associate’d with each coefficieuit, sc'e 
e*(|s. (7), (A 2 b), (A3b), and (A 5 b)) with n*sj)(*ct to 
the given gravitational coefficient. Since the* rates 
of change in perilune* altitude vary line*arly witli the 
gravitational coefficients, the sensitivity co(*fficie*nts 
are independent of the values of the coe'fficie'nts in 
a gravitational model. As a result, figures 26 to 28 
can be interpreted in two ways: Tlu'se figures re*pre*- 
sent a ‘'normalized” rate of change in perilune* alti- 
tude, independent of the value of the gravitational co- 
efficients, or they serve as a means of illustrating the* 
sensitivity of the perilune altitude decay rate with 
respect to the uncertainty in the value of the^ partic- 
ular co(*ffici(*nt. The contours of tlic'se figures are* in 
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Argument of perilune, co, deg Argument of perilune, co, deg Argument of perilune, co, deg 




(h) 



Fif>,ur(’ 28. Contour plot of s(*nsitivity for f ’o in units of 
kin/day x 1()‘^ for hp -- 100 kin. 


units of km/(lay x 10 ^ Plots of the sensitivity ro- 
('fiidents for a 300-km initial perilum' alt itiide orbit 
exhibited th(' saim' pattern as for tlu' 100 -kin cas(\ 
although tlu' values of the sensitivity coefficacmts \wvvo 
roughly only half as large. 

To evaluate the error in orbital lifcdinie predie- 
tions associated with the Ferrari inod('b the' variance 
in th(' pt'rilune d(H‘ay rate was eakailatc'd with th(' 
following forimila: 


'Orn 








'Or 


1 /2 




* \ 




wh(U'(' i'p is lh(' perilune altitude' rate' of change for 
the' corre'sponding coe'fiicieni (the' gravitational co- 
('fficients in the nuiiu'rator e)f ('epiation (0) appe'ai' as 
subsea ijits e)f Vp. not as value's of the coe'fhcients t lu'in- 
se'lve's) and <r is the' corn'spe)neiing stanelarel eleviat ion 
of ('ach gravitational ce)effici('nt value. In using this 
forniula, no cx)rrelat ion in the uncertaint ies of t he' val- 
ue's of the gi'avitational coefficie'uts is assiime'd. How- 
ever, the anel J 5 coefficie'uts are' highly e'orre'lat e'el. 
Since the'se exiefhcients affect the' orbit al e'leine'nts in a 
similar manner (table 3), se'parating the' e'fleets that 
are inelividually cemtributed by the'se' coe'flieaents is 
difficult. By ne)t accounting for {‘orre'lation eth'e'ls, 
e'ejuatiejii (9) ewe're^stimatt's the erreu s assexiate'd wit h 
t he' uncc'rtainties in the value's e)f ./;p . 75 , and (if 
the' anel J 5 terms are e lorn in ant. the erre>rs may be 
s i gn i fi c an 1 1 y o ve res t i mat ed ) . 

Fe^r the' Ferrari mejdel (ref. 23), rr/,^ — 1.8 x 10 
— 2.0 X 10“'b — 1.9 X 10“^b Figure's 29(a) 

te) 29(c) she)w the variance in perilune decay rate' in 
units e)f km/day for Q — 0°, 45 °, anel 90^. The simi- 
larity in tlmse plots indie:ates that the' unex'rtainty in 
the value of the C 31 ('oeffiedent c‘e)ntribute's little' to 
the ove'rall errea* in orbital lifetime pre'elictions. The 
similarity between the variance plots anel t he se'nsi- 
tivity coefficient plot for Jr> (hg* ^’i") indicate's that 
the' e'lTor in orbital decay rate's for the Fe'rrari model 
is eliie' aluie)st e'litirely to tlie uncertainty in the value 
e)f the Jr> coc'fficient. In fact, the standard de'viat ie)u 
assoeuate'd with J 5 is an order of magnituele' large'r 
than the? staiidaref deviations for J 3 and 7^1- 

As with the figures of the perilune' altitude' eiecay 
rates (fig. 21 ). the variance plots conve'y e)iily instan- 
tane'e)us information. Te> asse\ss the uncertainty in 
the pe?rilune altitude at the end of an orbit's lifetime. 
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Argument of perilune, co, deg Argument of perilune, co, deg Argument of perilune, (o, deg 




Inclination, i, deg 
( 1 )) n - \r,\ 



(c) ii = mr. 

Figure 29. Contour plot of Vf, variance in kin/day versus i 
ami u,' for an initial hp = 100 krn. 



Initial inclination, i, deg 


(a) Int('grated variance in kilonu*l('rs. 



Initial inclination, i, deg 


(1)) Lifetime in days. 

Figure 30. Orbital liL'tirm' cha.ra.ct('ristics for lOO-krn initial 

pcailuiK' altitude and il = (F. 

the variance must he iiitcgratc<l over tiriuy The com- 
puter program was modified to integrate the variance 
in a manner similar to the integration of the orbital 
elenumts (eq. (8)); this integration resulted in fig- 
ure 30(a). To interpret the figure, one must also reh^r 
to the plot of orbital lifetimes (fig. 30(b)). For ex- 
ample, a parking orbit with i = 40® and = 120® 
has a lifetime of more than 360 days (fig. 30(b)). 
However, this orbit has an uncertainty of approxi- 
mately ±40 km in the final perilune altitud(' as a rev 
suit of the inaccuracy in the gravitational coefficients 
(fig. 30(a)). 

The unce'rtainty in perilune altitude is less for 
orbits with short lifetimes, because the e^rrors in 
the decay rates have less time to accumulate. The 
largest uncertainties in perilune altitude occur for 
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orbits with lifetimes longer than a year. The un- 
certainties in the final values of perilune altitude 
range from 20 km to 300 km; this range demonstrates 
the inadequacy of current models for making long- 
term predictions of the orbital elements. A decrease 
in the uncertainty in perilune altitude^ valiums by an 
order of magnitude (which can only be achieved by 
decreasing the uncertainty in the vahu^s of the grav- 
itational coefficients themselves) is required to pro- 
vide mission planners with acceptable orbital lifetime 
predictions over the time periods addressed in this 
analysis. 

Applications of Analysis 

Although the uncertainty of the coefficients in the 
Ferrari model creates difficulty in gem'rating long- 
term orbital predictions, this analysis illustrates a 
technique to address effects of the nonspherical mass 
distribution of the Moon on a low-altitud(^ near- 
circular orbit. The usefulness of introducing a sim- 
plified gravitational model wIkui performing orbital 
lifetime studies at the preliminary mission dt^sign 
level has also been demonstrated. The Ferrari grav- 
itational model was specifically chosen for the pur- 
pose of illustrating techniques for determining or- 
bital lifetimes to enabk' quantitative results to be 
generated. However, several criteria can l)e applied 
from this analysis to assess the qualitative^ (ffiV'cts of 
other models on orbital lifetinu's. Tlu^sc' critenaa in- 
dicate the important parameters involve'd in lifetime 
predictions. 

Investigation of the gravitational coeffidenit 
illustrates that off-diagonal terms may have a signif- 
icant ('ffect on periluiu^ altitude decay rates. How- 
ever, these effects are periodic, with a 27.3-day cy- 
cle, and have little inffuence' on long-term lifetime 
predictions. Therefore, the consequences of lunar 
gravity fields on long-term orbital predictions are 
dependent mainly on the odd zonal harmonics, l)e- 
caiise they have tlu' greatest effect on th(^ eccentricity 
rates. The odd zonal coefficients model the geomet- 
rical asymmetry between the northern and southern 
hemispheres. Although both hemisf)heres contain the 
same amount of mass, the distribution of the mass 
within the two hemispheres is not the sanu'. (See 
ref. 22.) This mass distribution results in a long- 
period variation in the eccxmtricity. Thesc^ variations 
are proportional to {R/ay\ where n is th(' degree 
of the zonal term. (See ref. 4.) For high-altitude 
orbits, the effects due to higher order terms ar(' at- 
tenuated. However, for low-altitude orbits (such as 
the ones considered in this analysis), R/a converges 
slowly and enables higher order terms to contribute 
significantly. The figures of the sensitivity coefficients 


(figs. 26 to 28) can be used to calculate the rate of de- 
cay of perilune altitude for any specified value of the 
gravitational coefficients. Although the rates must be 
intc'grated to determine orbital lifetimes, a compar- 
ison of the rates indicates which terms in the grav- 
itational model provide significant contributions to 
altitude decay and shows the effect of decay rates as 
a function of inclination. 

For the Ferrari model, the J5 coefficient con- 
tributed the primary influence in orbital lifetimes 
(and the primary source of error); this influence 
indicates that long lifetimes exist near inclinations 
of 40.09° and 73.43°, where the eccentricity rates 
for J5 are zero for near-circular orbits (eq. (A3b)). 
For other gravitational models, more than one gravi- 
tational coefficient may strongly influence orbital life- 
times. The eccentricity rate for J3 is zero for an 
inclination of 63.43°. These results indicate that 
if J3 and J5 are the dominant odd zonal terms of 
the gravitational field and are of the same sign, long- 
lifetime orbits will exist somewhere within a 10° in- 
clination band (between 63.43° and 73.43°, figs. 22 
and 23). However, if the sign of J5 is negative (the 
actual sign of J5 is still unknown), a long-lifetime 
band of orbits will exist somewhere betwexm incli- 
nation values of 40.09° and 63.43°. These types of 
low- altitude, near-circular orbits are idc'al for mission 
planners that want low rates of decay in perilune al- 
titude (corresponding to a low AV budget for alti- 
tude boost maneuvers) and high inclination. (High 
inclination orbits are ideal for obtaining maximum 
coverage for a mapping mission or for providing ren- 
dezvous capability over a large range of landing-site 
latitudes for a manned mission.) 

A general observation can also be made about the 
lifetimes of polar orbits that are prediett'd by various 
spherical harmonic models. As a function of incli- 
nation, the ^ value that corresponds to each odd 
zonal harmonic will have a local minirnuin or max- 
imum (dependent on the sign of the coefficient) at 
i — 90° (fig. 23). Furthermore, since the odd zonal 
harmonics solely (to a first-order approximation) de- 
termine the long-period rate of change of the per- 
ilune alt i tilde, selection of polar parking orbits may 
result in large AV station- keeping requirements, de- 
pendent on the signs of the coefficients. If the signs 
of individual coefficients in the mod(4 are such that 
the ^ raters reinforce each other, the AV penalties 
might be substantial. Conversely, if the signs of the 
coefficients are such that the individual cont ributions 
interfere destructively, polar orbits may exist that 
have long lifetimes. This observation explains the 
discrepancy in the results shown in figure 7. For (ex- 
ample, the signs of J3 and J5 for the Ferrari 5x5 
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model are such t hat t heir occent ricity rates construc- 
tively inti'rFcn' and yield short polar orbit lifctiine 
])]■('( lift ions. Conversely, the' sign and magnitude of 
th(’ odd /onals for the Bills-Ferrari 16 x 16 model al- 
low the Jy tc'rm's large (Huitrihution to the eccen- 
tricity rat(^ to effectively cancel out the contribu- 
tions of the other terms. This supc^rposition yields 
long-lifet iuK' ])r('dictions for polar orbits. 

Coiichiding Remarks 

Antici])ation of future missions involving the 
pla<*('ment of spacecraft in circular low-altitud(' park- 
ing orbits about the Moon for long pcTiods of time 
motivat('d an investigation of the effects of the lu- 
nar gravitational field on an orbiting satellite. The 
formulation of a simplified five-coefficient lunar grav- 
itational modc'l was derived from Ferrari 's 5 x 5 
gravitat ional modt^l for the purpose of aiding mis- 
sion plamua’s in implementing nonspherical gravita- 
tional ('tfVn'ts of the Moon in prefiminary lunar mis- 
sion studic's. Tlu" ]>roi)os('d simplified gravitational 
mod('h consisting of the coefficients J?, ./f,, C22« 

and ^ vvas adefpiate in j)redicting orbital lifetimes 
of lunar satcdlites; the results generated from its 
use closi'ly corrc'sponded to thos(^ gtaicaated by the 
fVrrari nicxU^l. Furt.lu'rmor(\ tlu^ simplified model 
i-(’(1u(‘(h 1 tlu^ amount of c’omputer time required for 
liietinu' prt'diHions l)y approximat('ly two orders of 
magnitude': this reduction greatly c'lihaiK'ed the ca- 
[)ability to invt'stigate the orbital lifc'tinu's of lunar 
satellite's as a functie)n e>f the various initial orbital 
pai“ai!U'te'i‘s, 

He'siilts ge'iK'rateei with the' Fe'rrari simj)lified five- 
co('liie“i('nt gravitational nuKle'l ineli(‘at('el the e'xis- 
teauH' of se've'ral ine*linatie)n banels e)f short-lifetime 
and long-life'titue' e)rbits. These' banels refie'ct tlie be- 
havior of t lie' ])e'rilune' altituele' ele'e*ay rate' that is asso- 
e’iate'el with the' J5 e-e)effieae'iit . the' ele)minant term in 
the' ferrari mode'L as a fune’tie)ii of inclinatiein. Le)iig 
oriiital lifelime's were pe)ssil)le for lenv. miellatituele. 
and high-ine linat ieni e)rbits. althe)ugh pe)lar e^rbits had 
re'lative'ly slieirt liletime's (le'ss than 180 eiays) aex'orel- 
ing to this moelel. Of particular iiiterc'st is the pre- 
die tion oi’ a narreiw banel e)f e)rbits with ineTnatiems 


between appre)ximate'ly 60 ° anel 75^^' that yie'ld high 
lifetimes while prewieiing high-ine*lination e>rbits that 
are elesiraiile for various missie)ns. The life'times also 
ele'pemd on the initial value of argume'nt e)f perilune. 
especially at low inclinatiejiis. 

The C;u gravitational term ge'iic'rate'el elramatic 
localizeel effe^cts in the behavmr e>f eadutal lifeitime's, 
but did not contribute any sexailar effects to long- 
term predictions of lifetimes. The j)ur{dy localized 
effecds due to off-diagonal gravitational coefficients 
suggest that accurate long-term orbital j)r{'dictions 
can ])C performed by consideration of only the zonal 
harmonics. S])ecifically. consid('ration of the odd 
zonal terms may be sufficif'iit, as they hav(' the 
primary influence on orbital-altitudc' secular d('cay 
rates. This technique provides a simple nu'thod 
to generate orbital lifetime predictions with other 
models. 

Finally, the large uncertainty in the vafiu's oi“ the 
coefficients for the Ferrari model (or any other lunar 
gravity model) indicates that little confid('n(“(' can Ix' 
placed in the results tliat are generated. In partic- 
ular, the uncertainty in the value aiui evc'ii the sign 
of the J5 coefficit'iit (a rc'sult of corrc'laticjii between 
the J3 and the J5 term) contributes tlu' dominant 
error aSvSociated with tlu' F('rrari model. Currc'iit 
capabilities for long-term predictions of orbital life- 
tinu's for lunar satellites leave' much to Ix' dc'sin'd. 
Howc'ver, tlu' nu't hods j)r('sented in this analysis an' 
Ix'neficial for inc‘orporating tlx' Moon's nonsplu'rical 
gravitational ('ff('(‘ts on tlu' preliminary de'sign k'vx'l 
for future lunar mission planning with litth' addi- 
tional computational time ix'cjuin'd. Fur t lie r work 
needs to be performed in the determination of an ac- 
curate lunar gravity model, as cum*nt models {'ither 
give' inconsistent preelictions or pre'dictions with such 
large uncertainty value's that use'ful or me'aningful 
inte'rpre'tation of the' re'sults is difficult . 
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Tabl(‘ 1. Nuiiiorical Values (Uiiiioniializcd) of Coeffici(nils in Fc'rrari Model 
[Boldface' imnilH'rs denote values used in siniplilic'd mode'!] 
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-0.00789 

4 
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1 
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5 
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5 

3 

-0.00148 
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5 

4 

0.000598 

0.000456 

5 

5 
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Tabic 2. Numerical Value 

s (Unnormalized) of Coefficients in Bills-Ferrari Model 
[From ref. 32] 

Degree 

Order 

C'nm 

Sum 

2 

0 

-2.024 X 10^4 


2 

1 

-0.0904 X 10^^ 

0 

2 

2 

0.2226 X 10“ ‘ 

0.1936 X 10 

3 

0 

-0.889 X 10^"' 


3 

1 

0.2372 X 10^ ‘ 

0.7160 X lO"'* 

3 

2 

0.0483 X lO 'i 

0.1626 X 10"'' 

3 

3 

0.2212 X 10^''" 

-0.3415 X 10 * 

4 

0 

0.1173 X 10“^ 


4 

1 

-0.4573 X 10^^ 

0.1812 X 10"'' 

4 

2 

-0.1818 X lO-'’ 

-0.1512 X 10 

4 

3 

0.2868 X 10“'^ 

-0.8623 X 10"** 

1 

4 

-0.7396 X 10“’^ 

-0.1162 X 10 ^ 

5 

0 

-0.2388 X lO""’ 


5 

1 

-0.8272 X 10“'’’ 

-0.1310 X 10 ■'' 

5 

2 

0.6003 X 10“*" 

-0.3802 X 10"** 

5 

3 

-0.1287 X 10“^ 

0.1622 X 10 ** 

5 

4 

0.4360 X 10“*^ 

-0.5123 X 10 ' 

5 

5 

-0.1647 X 10-'^ 

0.2856 X 10"' 

f) 

0 

0.1774 X 10^* 


6 

1 

0.9127 X 10-^® 

-0.5508 X 10"' 

6 

2 

-0.5735 X 10-® 

-0.4491 X 10 ** 

6 

3 

-0.4830 X 10-^ 

-0.8541 X 10"' 

() 

4 i 

0.7759 X 10~^ 

-0.2914 X 10"** 

() 

5 

0.2663 X 10“*^ 

-0.6699 X 10"** 

G 

6 

0.1824 X 10“^ 

-0.2227 X 10"** 

7 

0 

0.2227 X 10^* 


7 

1 

0.1742 X lO-"’ 

0.1632 X 10 

7 

2 

-0.1753 X 10^** 

-0.4582 X 10"' 

7 

3 

-0.8736 X 10-*^ 

0.1479 X 10"' 

7 

4 

0.1890 X 10-^ 

0.7858 X 10"’* 

7 

5 

0.7750 X 10-*’ 

0.3291 X 10 ** 

7 

6 

-0.1374 X lO-'J 

0.4588 X 10"’* 

7 

7 

-0.2059 X 10-^ 

0.1702 X 10 ’* 

8 

0 

0.1493 X l()- ‘ 


8 

1 

-0.7009 X U)-*' 

0.1752 X lO"*’ 

8 

2 

0.1840 X 10“^ 

-0.4353 X 10"'’ 

8 

3 

-0.7987 X 10^*^ 

0.1688 X 10"'’ 

8 

4 

0.588 X 10"^ 

-0.1436 X 10"** 

8 

5 

-0.600 X lO'*-^ 

0.78 X lO"** 

8 

6 

-0.111 X 10"’* 

-0.104 X lO"'-* 

8 

7 

-0.153 X 10" *‘* 

-0.15 X 10"" 

8 

8 

0.23 X 10"*' 

-0.97 X 10"*'^ 







Appendix A 


First-Order Equations for Rates of Change of Orbital Elements 

First -or(l('r e(|uati(ms for rate's of change of the orbital eleiiieiits due to J2, </r>, C\\\ 

as follows: 
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(A5e) 


where n'pn^sents the ^ contribution for the C 31 term (eq. (A5d)). The ^ term represents the rate 

of change in the inertial noci(^ longitude of the orbit ft referenced to the prime meridian at epoch (the prime 
meridian is (hdiiu'd by th(' Moon- Earth direction at some specified time). The sclenographic node longitude of 
th(^ orbit Us is ndated to the inertial node longitude by the expression H.s = Uj - ivt. (See fig. AT) 

Idu'si' (^Illations (refs. 3 and 11 ) were in the computer program that was used to generate the plots for the 
three-dinu'iisional and contour plots. They were also used in the sensitivity studies. 


Polar view' 




Figure AT Differentiation between ^1, and Us- 
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Appendix B 

Verification of Data Generated by 
Program LUNLIFE 

The lunar lifetime program LUNLIFE^ was origi- 
nally intended to be the primary computational tool 
used in the present analysis. For this reason, much 
effort was invested in the verification of the data gen- 
erated by this program. Although the development 
of the simplified gravitational model outlined in ap- 
pendix C lessened the need for the use of the program 
LUNLIFE in this analysis, this program was essen- 
tial for the formulation and verification of the sim- 
plified model, LUNLIFE was also used to quantify 
the effects due to solar radiation pressure (including 
shadowing effects) and Earth-Sun gravity perturba- 
tions (utilizing an analytical ephemeris algorithm). 
A direct descendant of software used in Viking mis- 
sion studies, LUNLIFE numerically intc'grates the 
Lagrange planetary equations of motion (eqs, (4)) 
with a fixed-interval fourth-order predictor-corrector 
algorithm, and addresses perturbations through the 
use of disturbing functions developed by Kaula. (See 
refs. 40 and 41.) 

LUNLIFE was initially tested by comparing it to 
results contained in a Boeing report on a lunar model 
for Apollo. (Sec ref. 3.) The model in tlic^ Boeing 
report was ideal for comparison since it consisted of 
only four gravitational coefficients. A list of first- 
order analytical equations in the appendix of the 
Boeing report also allowed for testing of the effec'ts 
of individual gravitational coefficients on each of the 
orbital elements. Results generated by the Goddard 
Space Flight Center in the 1970's for a lunar lifetime 
study (ref. 4) were also duplicated with LUNLIFE, 


In this study, three 3x3 lunar gravity models were 
used to specifically address satellite lifetimes in near 
j)olar orbits. 

Additionally, LUNLIFE was compared in fig- 
ure B1 with the computer programs ASAP (ref. 39) 
and LUNNODE^ These two programs differ from 
LUNLIFE in that they integrate the equations of 
motion in Cartesian coordinates. As a result, much 
smaller time steps are used, and the programs in- 
clude short-period effects rather than averaging them 
out. The comparisons were made with an 8x8 
truncat(^d version of the Bills-Ferrari 16 x 16 grav- 
itational model. (See ref. 32.) The threes j)rograms 
generated near-identical data of perilune altitude ver- 
sus time. The various sources have validated the 
data that was generated by LUNLIFE and liave ver- 
ified the results that were generated by the medhod 
introduced in appendix C. 



Fif:;ure Bl. Verification of rc'sults gcncrate'd })y program 
LUM4FK. 


^ Developed by Flight Mechanics k Control under NASA 
contract NAS 1-1 8267 in February 1989. 
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Appendix C 

Coiiiparisoii of Simplified Gravity Model 
With Ferrari 5x5 Model 

llio fonmilation of a simplified fivc-coefficicnt lu- 
nar gra\dtati()iml model was derived from Ferrari’s 
5 X 5 gravitational model for pnri>ose of aiding 
mission planners in implementing nonsphcrical grav- 
itational effects of the Moon in preliminary lunar mis- 
sion studii's. The proposes I simplified gravitational 
niod(d consisted of the (*o(diici(aits J 2 . T5, C 22 ' 

and t'n- Several cases wine analy/(Hl with various 
initial orbital [)arameters to examine the changes in 
ptailnne altitude, longitnd(' of ascending node, and 
argiirrient of perilune between the Ferrari 5x5 model 
and the simplified five-coeffi(‘i(ait gravity model. 

A wide range of cases were analyzed by varying 
the initial inclination, longitude of ascending node, 
and argument of perilune. The initial mean anom- 
aly was set to z(a*o, the initial eccentricity w^is fixed 
at 0.05 (for ('onsideration of mvir-circular orbits), and 
tlu' scmiimajor axis value of 1935.79 km was selected 
to yield an initial perilune altitude of 100 km. Ini- 
tial inclination values of 45^ 90°, 120°, 150°, 
and 179° were selected to simulate both posigrade 
and i'('t rogradt^ orbits. In an (dfort to cover different 
([iiad rants, since cos a; and sin lc-’ appear in the ana- 
lytical ('xprc'ssions derivc'd from tlu' Lagrange piane- 
tary (‘(^nations, initial values of argument of pcrilinu' 
s(4ected were 0°, 135°, and 225°. Initial values of 
longitude" of asc’cnding node" of 0°, 135°, and 225° 
were" cliose'U fe>r the" same" re*ase)ns. The" pre"sene*e‘ e>l 
se"e“t()rial and te'sseral harme>nie’ coefficients in the" lu- 
nar nie)eiels intrexlucers a le)ngitudinal elepenelence in 
t lu" gravitat ional pe)te"ntial. This de*penelence" make*s 
it iKH'e'ssary te) utilize se"V("ral elilferent initial value's 
fe)r longit ude" e>f asce"nding nexiec If only ze)iial har- 
me)nie‘ e*e)e'ffi('ients we"re ineTide'd in the moelels, re"- 
sult.s we>ulei be" inelet)enelt"nt e)f initial values of longi- 
t iiele" e)f asce'uding nexle. Exhausting all ce)mbinatie)ns 
of the" values. 54 elifferent initial e‘e)nelitie)us were an- 
alyze"d. EaeT exise" was evaluatt'd for beTh the Ferrari 
5 X 5 !!K)elel and the" simplifie"el five-coefficient lunar 
gravity me)de"l. 

A single simplifie"el gravitat ie)nal model was pre- 
sunu'd te) he" insufficient fe>r ae‘e‘urate" e)rbital prexlie- 
tie)iis ove"r all ine'linations, l)("C'ause the magnitude" of 
the" ("ffe'e'ts e’e)iitril)uteel by each inelividual gravita- 
tie)iial co("ffieie‘nt is a functie)n erf iiudination. How- 
e"ve"r. the" simple me)del elid make accurate predic- 
tions over all latitudes that were te"sted. The r("sults 
of t he" Ih'rrari 5 x 5 model and the simplified mode"! 
we'H" analyze'el by e*e)m])aring their orbital lifetime's 


if they were less than 180 elays. If the lifetimes ex- 
(‘e"e"ele‘d 180 elays (e‘alculations were terniinate'd after 
180 days), the" twe) models w("re analyzexl by ce)in- 
paring their minimum i)erihiri("-altitueie value's elur- 
ing this time perie)el. Lifetime"s were the shortest fe>r 
inclinations of 90° and 120°, as each trial that was 
examined resulteel in imj)act with the hinar surface" in 
k"ss than 180 days. Near-e"(|uate)rial e)rbits (1°, 179°) 
had the smallest decrease in perilune" altitude; each 
case ("xam ineel had lifetime's in ("X(*e"ss e)f 180 days. 
These results are listed in tal)le" Cl. 

The data generate'd with the" simplified five- 
coefficient model matched very closely the data gen- 
erated with the Ferrari 5 x 5 mexlel. Fe)r cas("s with 
lifetimes of less than 180 days, the" simplified inoelel 
usually predicted the lifetimes within 4 days e)f the 
Ferrari model. Fe)r case's with lifetimes longer than 
180 days, the simplifie'el moelel usually predictexl 
the^ minimum perilune altit uele within 10 km e;f the 
Ferrari model. There? is iie> particular inclinatie>n 
range where the two gravitatie)ual rne)ele"ls se"("in to be 
in better agreement. Also, the simplifiexl me)del eloes 
iu)t consistently e)verestiniate or unde're?stiinate the 
orbital lifetimes (*e)ini)are'el with the Ferrari nie)eiel. 

There were two cases for whie:h the difference in 
predictions of the orbital lifetimes was signifi(‘ant. In 
the first case {i — 45°, U — 135°, uj — 135°), the sim- 
plified model underestimated the lifetime by 16 days. 
(See fig. CT(a).) The simplified model slightly under- 
predicted the rat(" of decay of tlu" perilune altitude 
through 74 days. This underprediction allowed the 
satellite to temporarily avoid a collision and, with 
the assistance of the medium-period effects, enabled 
the orbit to survive for an additional 16 days. In tlu" 
second case (/ = 45°, — 0°, — 225°). tlu" Ih'rrari 

mod("l predicted an orbital lifetinu’ of 129 days; the 
simplified model ])redicted a lih'time in excess of 
180 days, with a minimum {jcrilune altitude of 5 km 
for the first 180 days. (S("e fig. (M(b).) The explana- 
tion for this discrepancy is that the simplified model 
again underpredic’tcxl the rate of perilune dc'cay, al- 
lowing the orbit to exist past 129 days, and, with t he 
assistance of the long-period effects, enabU'd t he p("r- 
ihme altitude to actually increase. These two cases 
do not hinder the usefulness of the sim])lified model, 
as both of the discrepancies can be eliminated by 
assigning a small mai'gin of ("rror (arouml 5 km) to 
the value of the perilune altitiuk" predicted by t he 
simplified model. 

Also, the longitude of ascending node and ar- 
gument of perilune were monitored in this analy- 
sis. These orbital (?lenieiit values (conj])uted with the 
simplified and Ferrari models) also comparecl well, 
although not cjuit(" as well as the periluiH"-altitude 
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Perilune altitude, km Perilune altitude, km 



(a) Initial conditions: / = 45"^, — 145". u; = 145°. 



(b) Initial conditions: i — 15°, il = 0°, lj = 225°. 

Figure Cl. C'oniparison of fx^rilmH' altitmh' for siniplifit'd and 
Ferrari niod(4s. 

predictions. The loiigitudt' of ascending node preces- 
sion was consistent hetwt'en the two models, since tlie 
J ‘2 term, which is the dominant term for determining 


nodal precession, is included in both models. Agrtn^- 
nient of the argnnuait of perilune was not always as 
promising, since ruvir-circular orbits art' addri^ssed 
in this analysis. For near-circular orbits (wlu're a 
unique uj is not dt'hruHi), large changes in the argu- 
ment of perilune ofttui take place, these changt's t(md 
to amplify any dilfenuices between the two models. 

In th(' present analysis, the slniplifi(‘d gravita- 
tional niod(d generated results similar t,o thos(^ of tlu' 
Ferrari model for initial altitude's of 100 km. Veri- 
fication of accuracy for otluT initial altitude's woulel 
increase' the' applicability eT this moelel. The sinqdi- 
henl me)elel i>ropejseei in this analysis might be fur- 
ther impreweel by dire'ctly inap])ing the' e)bse'rvational 
eiata to the five coefliede'iits. as e)[)pe)seel te) the' pre'se'ut 
inetheKl (a truncateei version of the Ferrari me)ele'l) e>f 
aelopting the' value's e)f the e:e)e'ffiedents elire'edly fre)iii 
the Ferrari moele'l. 

This analysis suggests that e)t.he'r siniplihe'el grav- 
itatie)iial nu)elels coidd be pre)pose'd te> simulate' e've'ii 
higher e)rel('r moelels for the purpose of gene'rating 
orbital lifetime prexlictie^ns witlu)ut sacrificing accu- 
rate results. He)wever, the' capability of the' pro- 
pe)se'el nu)elel te) eluplie’ate' results similar te) the' Fe'i rari 
nioeiel has not yet bc'cn verifieel for other initial pe'r- 
ilune altitueles, anel the pre)i)e).se'el me)de'l is not e'x- 
pe'cted to generate results similar to t hose' pre'elicte'el 
by other gravitational moele'ls. Sincx' e)ther moele'ls 
eliffe'r in the magnituele' (anel pe'rha|)s eve'ii sign) of 
their coefficients, a selection e)f e'oeffiede'iits e)the'r than 
the one's state'd in the pre)pe)S('el simplifie'd moeh'l pre'- 
sented herein may be' nu)re' appre)priate' fe>r simulat- 
ing the results e)f that particular moele'l. Alse), ine)el- 
els with a larger iiuml)e'r e)f gravitatie)nal terms may 
be more difficult te> simulate with a simplifie'el niode'l 
be'cause e)f the magnituele e)f the' ex)iitribut ion of the' 
higher oreier terms. 
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Table Cl. Comparison of Orbital Lifetimes 


(a) i - 1°, 45°. and 90° 


Initial orbital elements Ferrari 5x5 model Simplified model 


L 

doK 

n. 

deg 

^7 

deg 

Orbital 

lifetime/^ 

days 

Minimum 

altitude. 

km 

Orbital 

lifetime/^ 

days 

Minimum 

altitude, 

km 

1 

0 

0 


42 


54 

1 

0 

135 


95 


92 

1 

0 

225 


78 


87 

1 

135 

0 


89 


88 

1 

135 

135 


00 


78 

1 

135 

225 


39 


54 

1 

225 

0 


84 


88 

1 

225 

135 


45 


50 

1 

225 

225 


73 


72 

45 

0 

0 


22 


25 

45 

0 

135 

G2 


05 


45 

0 

225 

129 



5 

45 

135 

0 


16 


13 

45 

135 

135 

74 


90 


45 

135 

225 


10 


10 

45 

225 

0 1 


13 


3 

45 

225 

135 

78 


70 


45 

225 

225 


18 


10 

90 

0 

0 

47 


47 


90 

0 

135 

102 


99 


90 

0 

225 

144 


145 


90 

135 

! 0 

45 


48 


90 

135 

1 135 

97 


97 


90 

135 

225 

139 


141 


90 

225 

0 

52 


48 


90 

225 

135 

105 


101 


90 

225 

225 

147 


143 



Orbitjil lifetime if < 180 davs. 




Table Cl. Concluded 


(b) i - 120°, 150°, and 179° 


Initial orbital elements Ferrari 5 x 5 model Sim{)lifi('d model 





Orbital 

Minimum 

Orbital 

Minimum 

i. 

ii. 

UJ ^ 

lifetime/^ 

altitud(\ 

lifetime," 

altitude. 

deg 

deg 

deg 

days 

km 

days 

km 

120 

0 

0 

167 


167 


120 

0 

135 

77 


77 


120 

0 

225 

77 


78 


120 

135 

0 

149 


(i>) 


120 

135 

135 

59 


59 


120 

135 

225 

60 


59 


120 

225 

0 

154 


148 


120 

225 

135 

44 


44 


120 

225 

225 

44 


46 


150 

0 

0 


5 


10 

150 

0 

135 


77 


75 

150 

0 

225 

136 


133 


150 

135 

0 


11 


4 

150 

135 

135 


46 


61 

150 

135 

225 

114 


114 


150 

225 

0 


19 


15 

150 

225 

135 


80 


75 

150 

225 

225 

127 


126 


179 

0 

0 


95 


96 

179 

0 

135 


66 


59 

179 

0 

225 


51 


63 

179 

135 

0 


53 


64 

179 

135 

135 


99 


99 

179 

135 

225 


74 


74 

179 

225 

0 


61 


65 

179 

225 

135 


71 


79 

179 

225 

225 


95 


94 


'Orbital lifotiriic if <180 days. 

Siiigularity in niod(d causc'd (‘rroiu'ou.s prcHliction. 
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